Overview

Explore Alaska Lightning data. Data for download within the zipped files on the site.

Lightning data to be used by Xanthe to evaluate wildfire probability.

Notes: – Temporal period: 2013 to 2022 – Spatial domain: Ecoregions of AK boreal interior, boreal cordillera, taiga plains, and taiga shield. – Lightning : yes I think we want it as precise as possible. However, in the proposal we outlined a subsampling plan where we get this data for 100 blocks/ecoregion stratified by fire occurrence (1 block ~ 50,000 ha). I am open to a different approach but maybe we don’t need lightning across the entire spatial domain? – 30m = yes – coordinate system = Canada Albers Equal Area Conic projection (ESRI:102001)

Notes on Lightning Data – Data have been clipped to show only those strikes that fall within 250 km of Alaska, Yukon Territory, British Columbia, and Northwest Territories # Packages

library(tidyverse)
## ── Attaching packages ─────────────────────────────────────── tidyverse 1.3.1 ──
## ✔ ggplot2 3.4.2     ✔ purrr   1.0.1
## ✔ tibble  3.2.1     ✔ dplyr   1.1.2
## ✔ tidyr   1.3.0     ✔ stringr 1.5.0
## ✔ readr   2.1.2     ✔ forcats 0.5.1
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag()    masks stats::lag()
library(sf)
## Linking to GEOS 3.10.2, GDAL 3.4.2, PROJ 8.2.1; sf_use_s2() is TRUE
library(lubridate)
## 
## Attaching package: 'lubridate'
## The following objects are masked from 'package:base':
## 
##     date, intersect, setdiff, union

Data

eco = st_read("../data/spatialData/dire_ecoregs_caea/dire_ecoregs_caea.shp", "dire_ecoregs_caea")
## Reading layer `dire_ecoregs_caea' from data source 
##   `/Volumes/Woodwell/Woodwell/R-GitHub/FiSL-Data-Production/data/spatialData/dire_ecoregs_caea/dire_ecoregs_caea.shp' 
##   using driver `ESRI Shapefile'
## Simple feature collection with 4 features and 8 fields
## Geometry type: POLYGON
## Dimension:     XY
## Bounding box:  xmin: -2939454 ymin: 1836187 xmax: 108924.9 ymax: 4333895
## Projected CRS: Canada_Albers_Equal_Area_Conic
shp = st_read("../data/lightningData/CurrentYearLightning_SHP/TOA_STRIKES.shp", "TOA_STRIKES")
## Reading layer `TOA_STRIKES' from data source 
##   `/Volumes/Woodwell/Woodwell/R-GitHub/FiSL-Data-Production/data/lightningData/CurrentYearLightning_SHP/TOA_STRIKES.shp' 
##   using driver `ESRI Shapefile'
## Simple feature collection with 598 features and 15 fields
## Geometry type: POINT
## Dimension:     XY
## Bounding box:  xmin: -169.3093 ymin: 51.26926 xmax: -135.4176 ymax: 70.3536
## Geodetic CRS:  WGS 84
txt = read_csv("../data/lightningData/CurrentYearLightning_TXT/TOA_STRIKES.txt")
## Rows: 598 Columns: 16
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr  (5): STROKETYPE, UTCDATETIME, LOCALDATETIME, STRIKETIME, POLARITY
## dbl (11): OID_, NETWORKCODE, MILLISECONDS, LATITUDE, LONGITUDE, AMPLITUDE, G...
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
historical = read_csv("../data/lightningData/Historical_Lightning_as_txt/Historical_Lightning_2012_2023_TOA_AlaskaAlbersNAD83.txt")
## Rows: 1834567 Columns: 16
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr  (5): STROKETYPE, UTCDATETIME, LOCALDATETIME, STRIKETIME, POLARITY
## dbl (11): OID_, NETWORKCODE, MILLISECONDS, LATITUDE, LONGITUDE, AMPLITUDE, G...
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
usa = st_read("../data/lightningData/cb_2018_us_state_20m/cb_2018_us_state_20m.shp", "cb_2018_us_state_20m")
## Reading layer `cb_2018_us_state_20m' from data source 
##   `/Volumes/Woodwell/Woodwell/R-GitHub/FiSL-Data-Production/data/lightningData/cb_2018_us_state_20m/cb_2018_us_state_20m.shp' 
##   using driver `ESRI Shapefile'
## Simple feature collection with 52 features and 9 fields
## Geometry type: MULTIPOLYGON
## Dimension:     XY
## Bounding box:  xmin: -179.1743 ymin: 17.91377 xmax: 179.7739 ymax: 71.35256
## Geodetic CRS:  NAD83

World map

https://datacatalog.worldbank.org/dataset/world-bank-official-boundaries

wrld = st_read("../data/spatialData/world/WB_countries_Admin0_10m.shp", "WB_countries_Admin0_10m")
## Reading layer `WB_countries_Admin0_10m' from data source 
##   `/Volumes/Woodwell/Woodwell/R-GitHub/FiSL-Data-Production/data/spatialData/world/WB_countries_Admin0_10m.shp' 
##   using driver `ESRI Shapefile'
## Simple feature collection with 251 features and 52 fields
## Geometry type: MULTIPOLYGON
## Dimension:     XY
## Bounding box:  xmin: -180 ymin: -59.47275 xmax: 180 ymax: 83.6341
## Geodetic CRS:  WGS 84

Dataframes

shp
## Simple feature collection with 598 features and 15 fields
## Geometry type: POINT
## Dimension:     XY
## Bounding box:  xmin: -169.3093 ymin: 51.26926 xmax: -135.4176 ymax: 70.3536
## Geodetic CRS:  WGS 84
## First 10 features:
##       STROKETYPE NETWORKCOD UTCDATETIM LOCALDATET MILLISECON LATITUDE LONGITUDE
## 1  GROUND_STROKE          0 2024-01-02 2024-01-02        327 60.24884 -145.3146
## 2  GROUND_STROKE          0 2024-01-03 2024-01-03        806 57.21603 -154.7037
## 3  GROUND_STROKE          0 2024-01-04 2024-01-03        514 59.64610 -147.5017
## 4  GROUND_STROKE          0 2024-01-04 2024-01-04        148 59.19926 -149.8828
## 5  GROUND_STROKE          0 2024-01-04 2024-01-04        419 58.24838 -159.9815
## 6  GROUND_STROKE          0 2024-01-04 2024-01-04         21 60.89884 -144.3514
## 7  GROUND_STROKE          0 2024-01-04 2024-01-04        708 57.26956 -152.1585
## 8  GROUND_STROKE          0 2024-01-04 2024-01-04        213 57.64578 -152.3733
## 9  GROUND_STROKE          0 2024-01-05 2024-01-05        839 60.64806 -144.8472
## 10 GROUND_STROKE          0 2024-01-05 2024-01-05        678 60.62782 -145.1794
##    AMPLITUDE GDOP ERRSEMIMAJ ERRSEMIMIN ERRELIPSEA       STRIKETIME STRIKESEQN
## 1     -53900    2          2          1        -84 2024/01/02 01:04      14601
## 2     -39700    8          9          1         33 2024/01/03 05:11      14633
## 3     -16700    2          2          1         82 2024/01/03 23:28      14809
## 4    -102900    2          3          1         46 2024/01/04 00:43      14858
## 5     -17300    8          9          1         16 2024/01/04 01:37      14881
## 6     -66900    2          3          1        -87 2024/01/04 04:52      14988
## 7     -26600    3          3          1          6 2024/01/04 05:21      15004
## 8     -39600    2          2          1          4 2024/01/04 13:18      15108
## 9     -98900    2          3          1        -87 2024/01/05 00:45      15135
## 10    -14300    2          2          1        -86 2024/01/05 00:48      15136
##    POLARITY                   geometry
## 1  Negative POINT (-145.3146 60.24884)
## 2  Negative POINT (-154.7037 57.21603)
## 3  Negative  POINT (-147.5017 59.6461)
## 4  Negative POINT (-149.8828 59.19926)
## 5  Negative POINT (-159.9815 58.24838)
## 6  Negative POINT (-144.3514 60.89884)
## 7  Negative POINT (-152.1585 57.26956)
## 8  Negative POINT (-152.3733 57.64578)
## 9  Negative POINT (-144.8472 60.64806)
## 10 Negative POINT (-145.1794 60.62782)
txt
## # A tibble: 598 × 16
##     OID_ STROKETYPE  NETWORKCODE UTCDATETIME LOCALDATETIME MILLISECONDS LATITUDE
##    <dbl> <chr>             <dbl> <chr>       <chr>                <dbl>    <dbl>
##  1     1 GROUND_STR…           0 1/2/2024 9… 1/2/2024 1:0…          327     60.2
##  2     2 GROUND_STR…           0 1/3/2024 1… 1/3/2024 5:1…          806     57.2
##  3     3 GROUND_STR…           0 1/4/2024 7… 1/3/2024 23:…          514     59.6
##  4     4 GROUND_STR…           0 1/4/2024 8… 1/4/2024 0:4…          148     59.2
##  5     5 GROUND_STR…           0 1/4/2024 9… 1/4/2024 1:3…          419     58.2
##  6     6 GROUND_STR…           0 1/4/2024 1… 1/4/2024 4:5…           21     60.9
##  7     7 GROUND_STR…           0 1/4/2024 1… 1/4/2024 5:2…          708     57.3
##  8     8 GROUND_STR…           0 1/4/2024 2… 1/4/2024 13:…          213     57.6
##  9     9 GROUND_STR…           0 1/5/2024 8… 1/5/2024 0:4…          839     60.6
## 10    10 GROUND_STR…           0 1/5/2024 8… 1/5/2024 0:4…          678     60.6
## # ℹ 588 more rows
## # ℹ 9 more variables: LONGITUDE <dbl>, AMPLITUDE <dbl>, GDOP <dbl>,
## #   ERRSEMIMAJOR <dbl>, ERRSEMIMINOR <dbl>, ERRELIPSEANGLE <dbl>,
## #   STRIKETIME <chr>, STRIKESEQNUMBER <dbl>, POLARITY <chr>
historical
## # A tibble: 1,834,567 × 16
##     OID_ STROKETYPE  NETWORKCODE UTCDATETIME LOCALDATETIME MILLISECONDS LATITUDE
##    <dbl> <chr>             <dbl> <chr>       <chr>                <dbl>    <dbl>
##  1     1 GROUND_STR…           0 3/7/2012 0… 3/6/2012 22:…          965     59.7
##  2     2 GROUND_STR…           0 3/7/2012 0… 3/6/2012 23:…          918     55.3
##  3     3 GROUND_STR…           0 3/7/2012 1… 3/7/2012 02:…          462     55.2
##  4     4 GROUND_STR…           0 3/7/2012 1… 3/7/2012 04:…          264     55.2
##  5     5 GROUND_STR…           0 3/7/2012 1… 3/7/2012 06:…          213     57.2
##  6     6 GROUND_STR…           0 3/7/2012 1… 3/7/2012 07:…          757     55.0
##  7     7 GROUND_STR…           0 3/7/2012 1… 3/7/2012 09:…          653     55.5
##  8     8 GROUND_STR…           0 3/8/2012 0… 3/7/2012 15:…          961     55.3
##  9     9 GROUND_STR…           0 3/8/2012 1… 3/8/2012 02:…          635     60.1
## 10    10 GROUND_STR…           0 3/8/2012 1… 3/8/2012 02:…          907     55.9
## # ℹ 1,834,557 more rows
## # ℹ 9 more variables: LONGITUDE <dbl>, AMPLITUDE <dbl>, GDOP <dbl>,
## #   ERRSEMIMAJOR <dbl>, ERRSEMIMINOR <dbl>, ERRELIPSEANGLE <dbl>,
## #   STRIKETIME <chr>, STRIKESEQNUMBER <dbl>, POLARITY <chr>
usa
## Simple feature collection with 52 features and 9 fields
## Geometry type: MULTIPOLYGON
## Dimension:     XY
## Bounding box:  xmin: -179.1743 ymin: 17.91377 xmax: 179.7739 ymax: 71.35256
## Geodetic CRS:  NAD83
## First 10 features:
##    STATEFP  STATENS    AFFGEOID GEOID STUSPS         NAME LSAD        ALAND
## 1       24 01714934 0400000US24    24     MD     Maryland   00  25151100280
## 2       19 01779785 0400000US19    19     IA         Iowa   00 144661267977
## 3       10 01779781 0400000US10    10     DE     Delaware   00   5045925646
## 4       39 01085497 0400000US39    39     OH         Ohio   00 105828882568
## 5       42 01779798 0400000US42    42     PA Pennsylvania   00 115884442321
## 6       31 01779792 0400000US31    31     NE     Nebraska   00 198956658395
## 7       53 01779804 0400000US53    53     WA   Washington   00 172112588220
## 8       72 01779808 0400000US72    72     PR  Puerto Rico   00   8868896030
## 9       01 01779775 0400000US01    01     AL      Alabama   00 131174048583
## 10      05 00068085 0400000US05    05     AR     Arkansas   00 134768872727
##         AWATER                       geometry
## 1   6979966958 MULTIPOLYGON (((-76.04621 3...
## 2   1084180812 MULTIPOLYGON (((-96.62187 4...
## 3   1399985648 MULTIPOLYGON (((-75.77379 3...
## 4  10268850702 MULTIPOLYGON (((-82.86334 4...
## 5   3394589990 MULTIPOLYGON (((-80.51989 4...
## 6   1371829134 MULTIPOLYGON (((-104.0531 4...
## 7  12559278850 MULTIPOLYGON (((-123.2371 4...
## 8   4922382562 MULTIPOLYGON (((-65.34207 1...
## 9   4593327154 MULTIPOLYGON (((-88.46866 3...
## 10  2962859592 MULTIPOLYGON (((-94.61792 3...
na = wrld %>% filter(FORMAL_EN %in% c("Canada", "United States of America"))

Explore types of data

unique(historical$STROKETYPE)
## [1] "GROUND_STROKE" "CLOUD_STROKE"
unique(historical$POLARITY)
## [1] "Negative"       "Positive"       "Cloud To Cloud" NA

Get AK boundary

ak = usa %>% filter(STUSPS=="AK")

Split Date time column into two

( historical1 = historical %>% 
  separate(LOCALDATETIME, into = c("LocalDate", "LocalTime"), sep = " (?=[^ ]+$)", remove = FALSE) %>%
  mutate(LocalDate = mdy(LocalDate)) %>%
  dplyr::mutate(year = lubridate::year(LocalDate), 
                month = lubridate::month(LocalDate), 
                day = lubridate::day(LocalDate))
)
## # A tibble: 1,834,567 × 21
##     OID_ STROKETYPE   NETWORKCODE UTCDATETIME LOCALDATETIME LocalDate  LocalTime
##    <dbl> <chr>              <dbl> <chr>       <chr>         <date>     <chr>    
##  1     1 GROUND_STRO…           0 3/7/2012 0… 3/6/2012 22:… 2012-03-06 22:43:29 
##  2     2 GROUND_STRO…           0 3/7/2012 0… 3/6/2012 23:… 2012-03-06 23:39:24 
##  3     3 GROUND_STRO…           0 3/7/2012 1… 3/7/2012 02:… 2012-03-07 02:40:50 
##  4     4 GROUND_STRO…           0 3/7/2012 1… 3/7/2012 04:… 2012-03-07 04:50:09 
##  5     5 GROUND_STRO…           0 3/7/2012 1… 3/7/2012 06:… 2012-03-07 06:47:30 
##  6     6 GROUND_STRO…           0 3/7/2012 1… 3/7/2012 07:… 2012-03-07 07:36:31 
##  7     7 GROUND_STRO…           0 3/7/2012 1… 3/7/2012 09:… 2012-03-07 09:23:14 
##  8     8 GROUND_STRO…           0 3/8/2012 0… 3/7/2012 15:… 2012-03-07 15:18:05 
##  9     9 GROUND_STRO…           0 3/8/2012 1… 3/8/2012 02:… 2012-03-08 02:15:53 
## 10    10 GROUND_STRO…           0 3/8/2012 1… 3/8/2012 02:… 2012-03-08 02:29:14 
## # ℹ 1,834,557 more rows
## # ℹ 14 more variables: MILLISECONDS <dbl>, LATITUDE <dbl>, LONGITUDE <dbl>,
## #   AMPLITUDE <dbl>, GDOP <dbl>, ERRSEMIMAJOR <dbl>, ERRSEMIMINOR <dbl>,
## #   ERRELIPSEANGLE <dbl>, STRIKETIME <chr>, STRIKESEQNUMBER <dbl>,
## #   POLARITY <chr>, year <dbl>, month <dbl>, day <int>

Quick view Plots

plot(st_geometry(shp))

plot(st_geometry(ak))

Create points shapefile for Historic Lightning

Function to make points

df_to_sf <- function(x){
  st_as_sf(x, coords = c("LONGITUDE","LATITUDE"), crs = 4326, remove = FALSE)
}

Drop Lat Lon NA and convert to shapefile

( historical2 = historical1 %>% df_to_sf() )
## Simple feature collection with 1834567 features and 21 fields
## Geometry type: POINT
## Dimension:     XY
## Bounding box:  xmin: -179.9885 ymin: 46.11086 xmax: 179.971 ymax: 79.37176
## Geodetic CRS:  WGS 84
## # A tibble: 1,834,567 × 22
##     OID_ STROKETYPE   NETWORKCODE UTCDATETIME LOCALDATETIME LocalDate  LocalTime
##  * <dbl> <chr>              <dbl> <chr>       <chr>         <date>     <chr>    
##  1     1 GROUND_STRO…           0 3/7/2012 0… 3/6/2012 22:… 2012-03-06 22:43:29 
##  2     2 GROUND_STRO…           0 3/7/2012 0… 3/6/2012 23:… 2012-03-06 23:39:24 
##  3     3 GROUND_STRO…           0 3/7/2012 1… 3/7/2012 02:… 2012-03-07 02:40:50 
##  4     4 GROUND_STRO…           0 3/7/2012 1… 3/7/2012 04:… 2012-03-07 04:50:09 
##  5     5 GROUND_STRO…           0 3/7/2012 1… 3/7/2012 06:… 2012-03-07 06:47:30 
##  6     6 GROUND_STRO…           0 3/7/2012 1… 3/7/2012 07:… 2012-03-07 07:36:31 
##  7     7 GROUND_STRO…           0 3/7/2012 1… 3/7/2012 09:… 2012-03-07 09:23:14 
##  8     8 GROUND_STRO…           0 3/8/2012 0… 3/7/2012 15:… 2012-03-07 15:18:05 
##  9     9 GROUND_STRO…           0 3/8/2012 1… 3/8/2012 02:… 2012-03-08 02:15:53 
## 10    10 GROUND_STRO…           0 3/8/2012 1… 3/8/2012 02:… 2012-03-08 02:29:14 
## # ℹ 1,834,557 more rows
## # ℹ 15 more variables: MILLISECONDS <dbl>, LATITUDE <dbl>, LONGITUDE <dbl>,
## #   AMPLITUDE <dbl>, GDOP <dbl>, ERRSEMIMAJOR <dbl>, ERRSEMIMINOR <dbl>,
## #   ERRELIPSEANGLE <dbl>, STRIKETIME <chr>, STRIKESEQNUMBER <dbl>,
## #   POLARITY <chr>, year <dbl>, month <dbl>, day <int>, geometry <POINT [°]>

write to shapefile

Change projections

roads <- st_transform(roads, crs = 4326) EPSG:3338

ls24 = st_transform(shp, crs=3338)
ak1 = st_transform(ak, crs=3338)
na1 = st_transform(na, crs=3338)
historical3 = st_transform(historical2, crs=3338)

FIlter Historical by summer month and years

historical13 = historical3 %>% filter(year==2013) %>% filter(month %in% c(5, 6, 7, 8, 9))
historical14 = historical3 %>% filter(year==2014) %>% filter(month %in% c(5, 6, 7, 8, 9))
historical15 = historical3 %>% filter(year==2015) %>% filter(month %in% c(5, 6, 7, 8, 9))
historical16 = historical3 %>% filter(year==2016) %>% filter(month %in% c(5, 6, 7, 8, 9))
historical17 = historical3 %>% filter(year==2017) %>% filter(month %in% c(5, 6, 7, 8, 9))
historical18 = historical3 %>% filter(year==2018) %>% filter(month %in% c(5, 6, 7, 8, 9))
historical19 = historical3 %>% filter(year==2019) %>% filter(month %in% c(5, 6, 7, 8, 9))
historical20 = historical3 %>% filter(year==2020) %>% filter(month %in% c(5, 6, 7, 8, 9))
historical21 = historical3 %>% filter(year==2021) %>% filter(month %in% c(5, 6, 7, 8, 9))
historical22 = historical3 %>% filter(year==2022) %>% filter(month %in% c(5, 6, 7, 8, 9))

Maps

Map area zoom

#### zoom to kamchatka, different projection ####

zoom_to <- c(-150.46, 64.93)  # ~ center of Kamchatka
zoom_level <- 3.5

# Lambert azimuthal equal-area projection around center of interest
#target_crs <- polarProj#sprintf('+proj=laea +lon_0=%f +lat_0=%f', zoom_to[1], zoom_to[2])


C <- 40075016.686   # ~ circumference of Earth in meters
x_span <- C / 2^(zoom_level+.1)
y_span <- C / 2^(zoom_level+.5)   # also sets aspect ratio

zoom_to_xy <- st_transform(st_sfc(st_point(zoom_to), crs = 4326), crs = 3338)
zoom_to_xy
## Geometry set for 1 feature 
## Geometry type: POINT
## Dimension:     XY
## Bounding box:  xmin: 167338.1 ymin: 1667711 xmax: 167338.1 ymax: 1667711
## Projected CRS: NAD83 / Alaska Albers
## POINT (167338.1 1667711)
disp_window <- st_sfc(st_point(st_coordinates(zoom_to_xy - c(x_span / 2, y_span / 2))),
                      st_point(st_coordinates(zoom_to_xy + c(x_span / 2, y_span / 2))),
                      crs = 3338)

Map Current 2024

Current Lightning is available as a shapefile for download.

Map 2013

Lightning May-September n=111,570

Map 2014

Lightning May-September n=116,352

Map 2015

Lightning May-September n=220,140

Map 2016

Lightning May-September n=232,625

Map 2017

Lightning May-September n=223,934

Map 2018

Lightning May-September n=76,987

Map 2019

Lightning May-September n=184,947

Map 2020

Lightning May-September n=108,943

Map 2021

Lightning May-September n=85,028

Map 2022

Lightning May-September n=169,551